Discrete Approximations for Singularly Perturbed Boundary Value Peoblems with Parabolic Layers, I∗1)

In this series of three papers we study singularly perturbed (SP) boundary value problems for equations of elliptic and parabolic type. For small values of the perturbation parameter parabolic boundary and interior layers appear in these problems. If classical discretisation methods are used, the solution of the finite difference scheme and the approximation of the diffusive flux do not converge uniformly with respect to this parameter. Using the method of special, adapted grids, we can construct difference schemes that allow approximation of the solution and the normalised diffusive flux uniformly with respect to the small parameter. We also consider singularly perturbed boundary value problems for convectiondiffusion equations. Also for these problems we construct special finite difference schemes, the solution of which converges ε-uniformly. We study what problems appear, when classical schemes are used for the approximation of the spatial derivatives. We compare the results with those obtained by the adapted approach. Results of numerical experiments are discussed. In the three papers we first give an introduction on the general problem, and then we consider respectively (i) Problems for SP parabolic equations, for which the solution and the normalised diffusive fluxes are required; (ii) Problems for SP elliptic equations with boundary conditions of Dirichlet, Neumann and Robin type; (iii) Problems for SP parabolic equation with discontinuous boundary conditions. General Introduction Consider a substance (or admixture) in a solution with a flux satisfying Fick’s law, and with distribution given by a diffusion equation. Let the initial concentration of ∗ Received March 4, 1995. 1) This work was supported by NWO through grant IB 07-30-012. 72 P.A. FARRELL, P.W. HEMKER AND G.I. SHISHKIN the admixture in the material as well as the concentration of the admixture on the boundary of the body be known. It is required to find the distribution of admixture in the material at any given time and also the quantity of admixture (that is the diffusive flux) emitted from the boundaries into the exterior environment. Such problems are of interest in environmental sciences in determining the pollution entering the environment from manufactured sources, such as houses, factories and vehicles, and from industrial and agricultural waste disposal sites, and also in chemical kinetics where the chemical reactions are described by reaction-diffusion equations. In considering such problems, it is important to note that the diffusion Fourier number, which is given by the diffusion coefficient of the admixture in materials, can be sufficiently small that large variations of concentration occur along the depth of the material. For small values of the Fourier number, diffusion boundary layers appear. Therefore these problems exhibit a singularly perturbed character. The mathematical formulation of such problems have a perturbation parameter which is a small coefficient (the diffusion Fourier number) multiplying the highest derivatives of the differential equation. Even in the case where only the approximate solution of the singularly perturbed boundary value problem is required, classical numerical methods, such as finite difference schemes and finite element methods[15, 16, 17] exhibit unsatisfactory behaviour. This arises because the accuracy of the approximate solution depends inversely on the perturbation parameter value and thus it deteriorates as the parameter decreases. In [18] it was shown that the use of classical numerical methods does not give approximate solutions with acceptable accuracy even for very fine grids. Thus, even the use of computers with extremely large capacity will not guarantee acceptable accuracy in the answer. To be more precise, it can be shown that the error in the approximate solution on any arbitrarily fine grid is greater than some positive number (independent of the number of grid nodes), for a sufficiently small value of the perturbation parameter (the diffusion Fourier number). For some applied problems such solution accuracy can be satisfactory. However even in these cases dissatisfaction can be caused by the lack of a guarantee than the use of a finer grid will increase the accuracy of the approximation. More serious problems occur when an accurate approximation of the spatial derivatives of the solution is also required. For example, in order to determine the quantity of admixture which enters the environment per unit of time, it is necessary to compute the gradient of the concentration of the substance along the normal to the surface of the material. When classical finite difference schemes are used it can be expected that errors in the computed diffusive flux will be much larger than those of the computed concentration. Such errors in evaluating fluxes can be often of unacceptable magnitude. Similar difficulties appear also in problems of heat exchange in cases where the heat Fourier number can take any arbitrary small value. One often requires an accurate approximation of the thermal flux on a boundary of the body. This series of papers is devoted to the construction of numerical approximations, Discrete Approximations for Singularly Perturbed Boundary Value Problems...... 73 using finite difference schemes, of singularly perturbed boundary value problems for elliptic and parabolic equations. The simplest example of problems of such type in a one-dimensional case is the problem of a stationary diffusion process with a reacting substance : ε d2 dx2 u(x)− c(x)u(x) = f(x), x ∈ D, u(x) = φ(x), x ∈ Γ , for c(x) ≥ c0 > 0, x ∈ D. Here D = ( 0, 1 ), Γ = D \D is a boundary of the domain D and the parameter ε can take any value in the interval ( 0,1 ]. The parameter ε characterises the diffusion coefficient of the substance and the function c(x) characterises the intensity of decay of the diffusion matter. When the parameter tends to zero, diffusion boundary layers appear in a neighbourhood of the boundary. In the case of regular boundary value problems the error in the approximate solution produced with the use of grid methods, is a function of the smoothness of the solution and of the distribution of the nodes of the grids used. However, the application of classical grid methods for such singularly perturbed boundary value problems leads to loss of accuracy for the approximate solution when the parameter value is small (see, for example, [12, 18] and results in the next section). The following question therefore arise: how to construct and to analyse special numerical methods for solving singularly perturbed boundary value problems, the approximate solution of which converges uniformly with respect to the parameter ε (or, in short, ε-uniformly). The error of the approximate solution obtained by such methods, should be independent of the parameter value and defined only by the number of nodes of the grid used. Detailed analytic investigations of such special numerical methods dates back to the end of 1960s (see, for example, [3, 12]). These first strong results for problems with boundary layers belong to two different approaches which are used for construction of special numerical methods: (a) fitted methods[12] on meshes with arbitrary distribution of nodes (for example, on a uniform mesh) the coefficients of difference equations (difference approximations) are chosen (fitted) to ensure parameter-uniform accuracy of the approximate solution; or (b) methods on special condensing grids (or adaptive meshes)[3]. Those methods use the standard classical difference equations but the nodes of the mesh are redistributed (or adapted, or condensed in the boundary layer) such that parameter-uniform convergence is achieved. Special, fitted schemes (that is the first approach) are attractive, since they allow the use of meshes with an arbitrary distribution of nodes, e.g. uniform grids (see, for example, [1, 2, 4, 6, 12] ). Using the second approach, adapted meshes with classical finite difference approximations, parameter-uniformly convergent schemes were also constructed for a series of boundary value problems (see, for example, [23] and references therein). For some boundary value problems parameter-uniformly convergent schemes 74 P.A. FARRELL, P.W. HEMKER AND G.I. SHISHKIN were constructed using either the first or the second approach for the same problem (see, for example [6, 20]), or using both approaches together for the same problem (see, for example [18, 19], where different approaches were used in different coordinate directions). In [14] both approaches were used at the same time (a fitted scheme on a grids with condensing nodes in the boundary layer). Thus there is a large variety of special approaches tailored to individual boundary value problems in the literature. In the case of singularly perturbed boundary value problems, for which accurate estimates of the diffusive fluxes are required, methods must be evolved which approximate both the solution and the normalised fluxes accurately. Investigations of such methods have been sparse in the literature (see, for example, [18]). In the present paper we consider singular perturbed elliptic and parabolic equations with parabolic boundary layers. For boundary value problems we construct special difference schemes, solutions of which converge ε-uniformly in an `∞norm. Also approximations of the normalised diffusive fluxes which converge ε-uniformly, are proposed. In the next section it is shown that the computed solution, for a singularly perturbed ordinary differential equation, which is found using a classical scheme does not converge ε-uniformly. We then consider the construction of special schemes which are ε-uniformly convergent. Grid approximations of solutions and diffusive fluxes for singularly perturbed parabolic equations are considered in this first paper. Approximations of elliptic equation with mixed boundary condition, which admit Dirichlet and Neumann conditions are studied in the second paper. To construct the special schemes in the first two papers methods based on special condensed grids are used. In the third and last paper we investigate singularly perturbed boundary value problems with discontinuous boundary conditions. In this case fitted methods are used. The improved special finite difference schemes which allow accurate approximation of both the solutions and the normalised diffusive fluxes for boundary value problems can be effectively applied for the solution and numerical analysis of applied problems with boundary and interior layers. The methods for construction of special schemes developed here can also be used to construct and investigate special schemes for more general singularly perturbed boundary value problems (see, for example, [7, 8, 23] ). The necessity to construct special schemes In order to demonstrate the problems which may appear in the numerical solution process, we consider the following simple example of a singularly perturbed ordinary differential equation for a boundary value problem: L(1.1)u(x) ≡ ε d2 dx2 u(x)− u(x) = −1, x ∈ D, (1.1a) u(0) = u(1) = 0; ε ∈ ( 0, 1 ], (1.1b) Discrete Approximations for Singularly Perturbed Boundary Value Problems...... 75 where D = ( 0, 1 ). For the solution of this problem we should like to use classical numerical methods, for example finite difference schemes2 . The standard scheme for the problem (1.1) is defined as follows. In the interval D, we introduce the grid Dh = ω1, (1.2) where ω1 is a uniform grid with a step-size h = 1/N , and N +1 is the number of nodes of the grid ω1. For the problem (1.1) we employ the classical difference scheme Λ(1.3)z(x) ≡ {εδx x − 1}z(x) = −1, x ∈ Dh, (1.3)